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ABSTRACT 

/ 

This  report  simmarizes  results  of  62  computations  performed  for 
and  included  in  the  proceedings  of  the  1981  Stanford  Conference 
on  Complex  Turbulent  Flows  (commonly  referred  to  as  the  l$8l 
Stanford  Olympics).  The  work  was  done  under  joint  sponsorship 
of  the  Air  Force  Office  of  Scientific  Research  (Contract  F49620-78- 
0-0024}  and  the  NASA  Lewis  Research  Center  (Contract  NAS3-228Q9). 

" " "  -- — - - 

-The  objective  of  this  study  has  been  to  use  a  single  set  of  equations 
modeling  turbulent  flow  phenomena,  with  no  adjustment  Of  closure 
coefficients  from  flow  to  flow*  to  predict  a  relatively  wide  range 
of  turbulent  flows*  In  so  doing  we  have  been  able  to  objectively 
assess  the  current  state  of  development  of  a  two-equation  model 
of  turbulence  and  to  establish  its  range  of  applicability,. 


Applications  include:  (a)  homogeneous  turbulent  flows;  (b)  incom¬ 
pressible  external  and  internal  flows;  and  (c)  compressible  external 
flows*  One  of  the  incompressible  cases  is  flow  past  a  backward 
facing  step  and  includes  boundary-layer  separation;  all  other  cases 
have  no  separation. 


The  model  employed  in  our  computations  predicts  flow  properties  in 
quite  close  agreement  with  experimental  data  for  the  constant- 
pressure  boundary  layer,  the  incompressible  mixing  layer  and  for 
flows  with  surface  mass  transfer.  Additionally,  predicted  effects 
of  Mach  number  and  surface  cooling  on  a  constant-pressure  boundary 
layer  are  close  to  measured  effects,, 

'-'T' 

For  flows  with  strong  adverse  pressure  gradient,  most  notably  the 
backward- facing  step,  the  model’s  predictions  differ  substantially 
from  corresponding  measurements.  As  an  almost  uniform  trend,  the 
model  initially  responds  to  an  adverse  pressure  gradient  similar  to 
what  has  been  measured  but,  upon  removal  of  the  gradient,  returns 
to  equilibrium  more  rapidly  than  measured. 
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1 .  INTRODUCTION 

The  objective  of  our  participation  in  the  1981  Stanford  Conference 
on  Complex  Turbulent  Flows  has  been  to  use  a  single  theory  to  com¬ 
pute  Is  wide  a  range  of  flows  included  in  the  Conference  as  possible. 
In  so  doing  it  his  been  our  hope  that  we  can  objectively  assess 
progress  made  to  date  in  developing  a  universally  applicable  en¬ 
gineering  model  of  turbulence.  To  accomplish  this  objective  we 
have  computed  20  flows  with  a  total  of  62  separate  computations. 

Table  1  summarizes  the  flows  computed,  including  the  Sponsor  for 
each  case. 


fable  1.  Summary  of  Flows  Computed 


All  b2  computations  used  the  two-equation  model  of  turbulence  de¬ 
vised  by  Wilcox  and  Rube sin"  with  some  minor  "fine  tuning"  of  the 
closure  coefficients.  Computational  tools  used  to  solve  the 
equations  of  motion  include:  (a)  a  fourth-order  accurate  Runge- 
Kutta  integration  scheme  for  the  8  homogeneous  turbulent  flows 
(Flows  0371,  0372,  037*  and  0376);  (b)  a  fully  elliptic  incompress 


program  named  EDDYNSI  for  the  backward-facing  step  (Flow  0421); 
arid  (c)  a  compressible/incOmpressible  boundary-layer  program 
named  EDDYBL  for  the  other  53-  computations . 

In  the  following  Sections,  we  present  in  detail  the  equations 
of  motion  and  boundary  conditions  employed.  We  then  give  a  brie 
description  of  the  numerical  tools  used  followed  by  a  discussion 
of  the  various  numerical  checks  made  during  the  course  of  the 
many  computations.  Finally,  we  summarize  results  obtained  and 
outline  possible  future  avenues  Of  research. 


2.  EQUATIONS  OF  MOTION 


This  section  first  presents  the  basic  equations  of  motion  used  in 
this  project.  Then  a  nonlinear  stress/strain-rate  constitutive 
relation  used  for  the  homogeneous  turbulent  flows  is  given.  Finally 
we  specify  Special  modifications  to  the  basic  model  needed  for 
flows  with  significant  system  rotation  and/or  streamline  curvature. 


2.1  CONSERVATION  EQUATIONS 


The  equations  of  motion  used  in  all  of  our  computations  are  those 

devised  by  Wilcox  and  Rubesin.  The  model  is  of  the  two-equation 

variety  in  which  the  Reynolds  stress  tensor  x. .  is  assumed  pro- 

portional  to  the  mean  strain  rate  tensor  S..  according  to 
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where  e  is  the  eddy  diffusivity,  e  is  the  turbulent  mixing  energy  * 
p  is  mass  density*  u^^  is  the  mean  velocity  vector,  x,  is  position 


vector  and  is  the  Kronecker  delta.  The  mean  equations  of  motion 


thus  are  witten  (for  steady  flow)  as  follows. 
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in  equations  (2-A)  p  is  mean  pressure,  h  is  mean  enthalpy,  u  is 
molecular  viscosity,  PrT  and  Prm  are  laminar  and  turbulent  Prand 
numbers,  cu  is  turbulent  dissipation  rate  and  S*  is  a  closure  co¬ 
efficient  which  will  be  defined  momentarily.  Before  introducing 


the  two  turbulence  model  equations  it  is  instructive  to  note  that 

the  mean  energy  equation  (4)  appears  to  find  the  conventional  work 

term  t,.,3u,,/3x.  replaced  by  3*pae.  This  is  not  an  ad  hoc  closure 
—  <3 

approximation,,  but  rather  a  closure  approximation  consistent  with 
those  made  below  in  the  turbulent  energy  equation.  The  correctness 
of  Equation  (k)  becomes  obvious  when  the  resultant  equation  for 
total  energy,  viz,  (h+%u.ui+e) ,  is  formed. 

To  complete  our  set  of  equations,  we  compute  the  eddy  aiffusivity 
in  terms  of  e  and  a  from: 

e  =  Y*e/ai  {5 

where  y*  is  a  closure  coefficient  given  below  in  Equations  (8). 

The  equations  governing  the  evolution  of  e  and  a  are: 
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where  l  is  turbulent  length  scale  defined  as  e^/a.  In  Equations 
(5-7)  there  are  several  closure  coefficients  whose  values  are  as 
given  in  Equations  ( 8 ) ^  , 

Y*  =  1  -  (1-a2)  exp  (-Re,-)  / 
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The  homogeneous  turbulent  flow  calculations  used  the  nonlinear  st 
strain-rate  constitutive  relation  devised  bv  Wilcox  and  Rube sin1 


which,  for  the  high  Re, 


incompressible  cases  considered  becomes: 


rpj 


0 


SeS, 


f^ij  +  9B*“(^+2S  3  )  (Simrimj  +  } 

mn  nm 


jm  mi' 


(9) 


where  H. ,  =  H  (3u./3x.  -  3u./3x. )  is  the  mean  rotation  tensor 
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2.3  SYSTEM  ROTATION  AND  STREAMLINE  CURVATURE 


For  rotating  homogeneous  turbulent  flow  (Case  0372)  commutations 

2 

include  the  Wilcox-Chambers  rotation  term.  This  term  is  added 
to  the  equation  for  e,  which  becomes 
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where  t  is  time  and  &  is  rotation  rate.  Finally,  the  transonic 

P 

airfoil  Cases  8621  and  8623  use  the  Wilcox-Chambers  streamline 
curvature  term  which  yields  the  following  modified  equation  for 
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where  R  is  surface  radius  of  curvature  (positive  for  convex, 
negative  for  concave). 


BOUNDARY  CONDITIONS 


For  all  but  the  3  homogeneous  cases  solid  boundaries  are  present 
and  surface  boundary  conditions  must  be  specified.  With  the  exceptio 
of  the  backward- facing  step,  all  computations  integrated  all  the 
way  to  the  surface,  y  =  0.  The  backward-facing  step  computation 
employed  "surface”  boundary  conditions  based  on  the  law  of  the  wall. 
This  Section  describes  the  surface  boundary  conditions  used  in 
each  case. 


i.l  INTEGRATION  TO  THE  SURFACE 


For  integration  all  the  way  through  the  sublayer  to  the  surface, 
y  =  0,  boundary  conditions  are  as  follows: 
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where  T  is  temperature  and  subscript  w  denotes  Surface.  For  the 
dissipation  rate,  the  quantity  S  is  a  universal  function  of  surface 
roughness  and  mass  injection  rate  defined  by  {Wilc6x-TraciJ) : 
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where  v  is  kinematic  viscosity. 


3.2  MATCHING  TO  THE  LAW  OF  THE  WALL 


For  the  backward-facing  step  we  used  boundary  conditions  consistent 
with  the  lav;  of  the  wall,  viz, 
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where  <  =  . 41  is  Harman’s  constant.  While  more-accurate  bcundarj 
conditions  (also  known  as  wall  functions)  are  available  for  our 
equations  (see  the  Appendix),  limited  time  and  funds  precluded 
their  use  for  the  Conference. 


4.  COMPUTATIONAL  TOOLS 
4-1  HOMOGENEOUS  TURBULENT  PLOWS 

Because  of  their  inherent  simplicity,  we  solved  the  equations  of 
motion  for  the  8  homogeneous  turbulent  flows  with  a  straightforward 
fourth-order  accurate  Runge-Kutta  integration  scheme  and,  for 
obvious  reasons,  the  program  requires  no  further  description. 

4.2  ATTACHED  AND  FREE  SHEAR  PLOWS 

The  lion's  share  of  our  computations  used  the  same  program,  namely 

our  two-dimensional/axisymmetric  compressibie/lncompressible 
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boundary-layer/shear-layer  program  known  as  EDDYBL  .  In  performing 
the  calculations  all  compressible  cases  were  done  on  a  UNIV-AC  1108 
and  all  incompressible  cases  on  a  TRC-80  Microcomputer.  The  latter 
cases  were  actually  done  with  a  version  of  EDDYBL  in  which  all  of 
the  compressibility  terms  were  eliminated.  The  program  is  a 
parabolic  marching  code  which  is  second-order  accurate  in  both 
streamwise  and  normal  directions  * 

4.3  BACKWARD  FACING  STEP 

The  backward-facing  step  case  was  done  with  an  incompressible, 
elliotic  program  known  as  EDDYNSI.  The  program  is  a  modified  ver- 
sion  of  the  TEACH-2E  Code  which  also  is  second-order  accurate  in 
streamwise  and  normal  directions. 

4.4  NUMERICAL  SENSITIVITY  STUDY 

We  performed  many  numerical  accuracy  tests  on  a  more  or  leSs  ran¬ 
dom  sampling  of  the  many  cases  we  computed.  In  general  we  tested  the 
effect  of  total  mesh  point  number,  location  of  mesh  point  nearest 
the  surfa  e  and  size  of  streamwise  steps  taken.  For  all  of  the 
boundary^-iayer  cases  we  found  80-100  mesh  points  normal  to  the 
surface  with  the  value  of  y+  for  the  point  newest  the  surface  less 
than  unity  to  be  quite  satisfactory.  Except  for  very  strong  adverse 
pressure  gradient  cases  there  is  virtually  no  loss  in  accuracy  in 
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taking  streamwise  steps  up  to  about  one  boundary-layer  thickness. 
In  some  of  our  compressible  boundary-layer  runs  we  used  as  many 
as  280  points  normal  to  the  surface  with  y  nearest  the  surface 
as  small  as  .09.  The  difference  in  computed  integral  properties 
over  a  100  point  calculation  was  never  found  to  be  more  than  2%. 

For  the  backward-facing  step  case  we  used  meshes  which  had  a 
total  of  196,  529 j  and  87 0  mesh  points.  The  total  number  of  mesh 
points  had  Very  little  effect  on  predicted  reattachment  length 
although  local  flow  properties  varied  substantially  with  the 
number  of  points  used. 
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5 •  RESULTS 


This  section  presents  a  case-by-case  description  of  results 
obtained  including  all  plots  submitted  to  the  1981  Stanford 
Olympics  Conference. 

5.1  HOMOGENEOUS  TURBULENT  PLOWS 

These  flows  have  no  sol  i  boundaries  and  diffusion  across  stream® 
lines  is  negligible.  n"-  s  the  equations  of  motion  simplify  to 
first-order  ordinary  dif*e  eotial  equations  (convective  terms 
are  replaced  by  time  rate  of  change  terms,  c.f.  Equation  10). 

The  equations  of  motion  are  trivially  integrable  using  a  standard 
Runge-Kutta  algorithm. 

Our  purpose  in  doing  these  flows  was  to  clearly  delineate  one  of 
the  bounds  on  the  applicability  range  of  a  two-equation  model  of 
turbulence,  The  assumption  of  an  algebraic  relation  between  the 
Reynolds  stress  tensor  and  the  mean  strain  rate  tensor  implies 
that  the  flow  has  achieved  an  "equilibrium"  state.  Even  using 
a  nonlinear  stress/strain-rate  constitutive  relation  (Equation  9) 
accounts  only  for  noneauipartition  of  energy;  departures  from 
equilibrium  still  are  ignored.  Results  for  the  four  homogeneous 
turbulent  flow  cases  follow. 

5.1.1  Homogeneous  Isotropic  Turbulence 

Figure  1  conucares  comouted  and  measured  turbulent  kinetic  energy, 

p 

q  for  decaying  isotropic  turbulence.  As  shown  computed  and  measured 
energies  differ  by  less  than  1%  of  scale.  This  close  agreement 
is  unsurprising  as  the  ratio  of  3  to  3*  has  been  selected  to  match 
measured  decay  rates  for  homogeneous  isotropic  turbulence, 

5.1.2  Homogeneous  Rotating  Turbulence 

Figures  2-5  compare  computed  and  measured  flow  properties  for 
three  different  rotation  rates,  vis,  Q=0,  20,  80  sec~^.  As  shown 
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in  Figure  2  (the  nonrotating  case),  we  agai'  predict  the  decay  of 

homogeneous  isotropic  turbulence  quite  accurately.  Sven  with 

-1  2 
20  sec  ,  predicted  and  measured  decay  of  q  are  quite  close.  How¬ 
ever,  at  the  highest  rotation  speed  we  actually  Dredict  an  eventual 
2 

increase  in  a  in  contrast  to  the  monotonic  decrease  measured. 

Figure  5  compares  computed  and  measured  ratios  of  /u 7 2  to  /v 1  * , 
Clearly  the  measured  partition  of  energy  differs  substantially 
from  that  predicted. 

5.1.3  Homogeneous  Plane  Strain 

—1 

Two  strain  rates  were  considered,  viz,  -3v/3y  =  3w/3z  =  9.44  sec 

(Townsend)  and  4.45  sec-1  (Tucker-Reynolds ) .  Figures  6-8  compare 

computed  and  measured  normal  Reynolds  stresses  for  the  higher 

strain  rate  while  Figures  9-11  correspond  to  the  lower  strain  rate. 

2 

As  shown,  for  both  cases,  w’  is  reasonably  close  while  predicted 
2-2 

u*  and  v’  are  about  50$  lower  than  measured. 

5.1.4  Homogeneous  Shear 

Two  shear  rates  were  considered,  viz,  3u/3y  =  12.9  sec-1  (Champagne, 
et  al)  and  48  sec-1  (Harris,  et  al).  Figures  12-15  compare  computed 
and  measured  Reynolds  stress  components  for  the  lower  shear  rate 
while  Figures  16-19  correspond  to  the  higher  shear  rate.  For  both 
cases  predicted  normal  and  shear  stress  components  are  much  lower 
than  measured. 

5-2  CONSTANT  PRESSURE  BOUNDARY  LAYER 

Our  next  round  of  applications  of  the  turbulence  model  is  to  flow 
over  a  flat  plate  at  both  incompressible  and  compressible  flow 
conditions.  In  all  computations,  computation  was  initiated  at  the 
leading  edge  of  the  plate  from  laminar  profiles.  The  equations 
of  motion  are  integrated  through  transition  up  to  Re^  =  10000 
for  all  of  the  compressible  cases. 
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5.2.1  Incompressible  Case 


Figures  20-22  compare  computed  and  measured  velocity  profile,  skin 
friction  and  shape  factor  for  an  incompressible  flat-plate 
boundary  layer  (FPBL).  The  velocity  profile  corresponds  to  a  plate- 
length  Reynolds  number  of  10.9  million.  As  shown,  differences  be¬ 
tween  computed  and  measured  flow  properties  are  well  within  engineer¬ 
ing  accuracy. 

5.2.2  Effect  of  Mach  Number 

Figure  23  compares  computed  effect  of  freestream  Mach  number  on 
an  adiabatic-wall  FPBL.  As  shown  for  Mach  number  ranging  from  0 
to  5  and  at  a  momentum^thickness  Reynolds  number,  Reg,  of  10000, 
the  model  equations  predict  skin  friction  approximately  3-5? 
lower  than  measured.  Figure  24  shows  the  predicted  recovery  factor 
as  a  function  of  Mach  number.  The  predicted  variation  is  well 
within  experimental  data  scatter. 

5.2.3  Effect  of  Surface  Cooling. 

Figure  25  compares  computed  and  measured  effects  of  surface  cooling 
on  a  Mach  5  FPBL.  The  adiabatic  wall  temperature  determined  from 
the  Mach  5  computation  of  Subsection  5.2.2  was  used  for  all  surface 
cooling  cases.  Again  computed  skin  friction  is  about  3-5?  lower 
than  measured. 

5.3  THE  MIXING  LAYER 

Perhaps  the  most  basic  of  all  free  shear  flows  is  the  mixing  layer. 
The  mixing  layer  is  the  next  of  our  applications.  In  this  subsection 
we  first  describe  our  results  for  the  incompressible  case,  Including 
effects  of  velocity  ratio.  Then  we  discuss  the  compressible  case. 

5.3.1  Incompressible  Case 

Application  of  the  Wilcox-Rubesin1  model  to  this  flow  found  the 
predicted  spreading  rate  to  be  .085  as  compared  to  the  measured 


12 


(consensus)  value  of  .115-  Further  investigation  showed  that  the 
spreading  rate  is  strongly  affected  by  the  closure  coefficients 
a  and  c*  (see  Equations  6-7).  Figure  26  shows  the  predicted 
effect  of  a  on  spreading  rate  (the  curve  was  constructed  with  a=b*). 

As  shown,  selecting  c  =  c*  =  2/3  yields  a  spreading  rate  of  .115* 

This  is  the  "fine  tuning"  of  the  VIlcox-Rubesin  model  alluded  to 
in  the  Introduction. 

Figure  27  compares  computed  and  measured  spreading  rate  as  a  function 
of  velocity  ratio  Up/u^  where  Ug  and  are  the  velocities  of  the 
mixing  streams.  As  shown,  predicted  spreading  rate  virtually 
duplicates  the  accepted  correlation  of  measured  values.  Figure  28 
compares  the  computed  velocity  profile  with  the  measurements  of 
Liepmann  and  Laufer. 

Finally,  Figure  29  compares  computed  and  measured  development  of  a 
mixing  layer  from  separation  to  a  distance  1800  momentum  thickness 
downstream.  The  predicted  asymptotic  spreading  rate  is,  as  expected, 
.115*  As  shown,  the  initial  spreading  rate  is  somewhat  higher  than 
measured  and  falls  a  bit  below  measured  values  farther  downstrea an 

5-3.2  Effect  of  Mach  Number 

To  assess  effects  of  compressibility,  we  next  compute  the  adiabatic 
mixing  layer,  viz,  the  mixing  of  a  supersonic  stream  with  a  stream  of 
the  same  fluid  at  rest  having  identical  total  temperatures.  In  order 
to  differentiate  effects  of  Mach  number  and  density  variation  we  first 
predicted  the  Mach  zero  spreading  rate  for  density  ratios  of  1/7  and 


7.  When  the  denser  fluid  is  at  rest  the  spreading  rate  is  reduced 
to  .111;  when  the  heavier  fluid  is  at  rest  the  spreading  rate  in¬ 
creases  slightly  to  .116*  Then,  varying  Mach  number  from  0  to  19, 
we  find  virtually  no  effect  whatever  on  spreading  rate  (Figure  30)* 

To  be  certain  no  numerical  errors  are  involved  we  reran  all  comp¬ 
ressible  cases  with  a  one-dimensional  implicit  time  marching  program 
to  solve  the  farfield  (self-similar)  equations  and  found  precisely 
the  same  result... no  variation  of  spreading  rate  with  Mach  number. 


Using  the  farfield  equations  we  also  included  (a)  transverse  pre- 
sure  gradient,  ( b 5  the  Saf f man-Wile ox^  compressibility  term  and 
(c)  the  Wilcox-Chambers  streamline  curvature  term.  None  of  these 
modifications  had  any  substantial  (greater  than  103)  effect  on 
predicted  spreading  rate. 

5.4  PLOWS  WITH  SURFACE  MASS  TRANSFER 

We  now  turn  to  more  advanced  applications,  the  first  of  which  are 
flows  with  surface  mass  transfer.  Only  incompressible  cases  'were 
done,  one  with  surface  mass  injection  and  the  other  six  computations 
with  Suction.  Results  follow. 

5.4.1  Effect  of  Blowing 

We  first  consider  a  constant  pressure  boundary  layer  with  uniform 
mass  injection,  F  *  v^/U^  =  .00375,  where  Ug  is  velocity  at  the 
boundary  layer  edge.  Figures  31-34  show,  respectively,  skin 
friction,  momentum  thickness,  displacement  thickness,  and  four 
velocity  profiles.  As  shown,  computed  and  measured  flow  properties 
are  quite  close.  The  largest  differences  are  in  momentum  and 
displacement  thickness  where  computed  differences  are  less  than  103. 

5.4.2  Effect  of  Suction 

Now  we  turn  to  suction  where  we  performed  a  total  of  six  computa¬ 
tions.  The  first  case  had  a  mild  adverse  pressure  gradient  and  a 
suction  rate  F  =  -.004,  Figures  35-38  compare  computed  and  measured 
skin  friction,  momentum  thickness  and  velocity  profiles.  Computed 
skin  friction  is  approximately  5%  higher  than  measured  while 
computed  and  measured  velocity  profiles  differ  by  less  than  7%,  Al¬ 
though  larger  differences  are  present  for  momentum  and  displacement 
thickness,  computed  and  measured  shape  factors  are  very  close. 

The  other  five  cases  are  all  at  the  same  freestream  flow  conditions 
and  have  zero  pressure  gradient.  Suction  rate  for  the  five  cases 
ranges  from  F  =  0  to  «*0l44.  Figure  39  compares  computed  and 
measured  velocity  profiles  for  the  unsucked  case  and  for  the  highest 
suction  rate.  As  shown,  differences  are  slight.  Figures  40-42 
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display  predicted  ana  measured  Reynolds  stresses  (as  with  the 
homogeneous  turbulent  flow  cases  we  used  Equation  9  to  compute  the 
normal  stresses).  At  the  three  highest  suction  rates  the  shear 
and  streamwise-normal  components  are  within  5%  of  their  measured 
counterparts.  The  lateral-normal  component  shows  larger  dis¬ 
crepancies. 

5*5  BOUNDARY  LAYERS  WITH  ADVERSE  PRESSURE  GRADIENT 

Our  attention  now  turns  to  effects  of  adverse  pressure  gradient, 
the.  long  standing  nemesis  of  turbulence  modelers.  Applications 
included  in  this  subsection  are  one  incompressible  computation 
and  ten  supersonic  computations. 

5*5*1  Incompressible  Case 

Our  incompressible  application  is  the  carefully  documented  flow 
of  Samuel  and  Joubert.  Figure  43  exhibits  skin  friction,  Figures 
44-45  show  shear  stress  profiles  and  Figure  46  displays  velocity 
profiles  *  As  a  general  observation,  the  computed  boundary  layer 
thickens  more  rapidly  than  measured  and,  rather  than  approaching 
separation,  tends  to  recover  more  rapidly  than  measured  as  the 
pressure  gradient  is  removed. 

5*5*2  Compressible  Gases 

Our  selection  of  compressible  cases  was  more  extensive  including 
a  round  of  nine  computations  at  Mach  2.2  with  three  different 
Reynolds  numbers  (Aeharaya,  et  al).  The  tenth  computation  was 
at  Mach  4  (Ewarts). 

Figures  47-49  compare  computed  and  measured  skin  friction  for  the 
nine  Mach  2.2  cases.  As  shown,  for  all  Reynolds  numbers,  the 
computed  skin  friction  begins  to  drop  at  about  the  same  location 
as  measured,  fails  to  drop  as  low  as  measured,  and  recovers  much 
more  rapidly  than  measured  as  the  adverse  gradient  is  removed. 
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Figures  50-51  show  shape  factor  for  two  of  the  cases.  As  illustrated 
the  measurements  indicate  rapid  variations  in  shape  factor  while 
the  predicted  shape  factors  vary  much  more  gradually.  Figures 
52-57  display  computed  and  measured  velocity,  turbulent  energy  and 
shear  stress  profiles.  As  in  the  incompressible  case,  the  numerical 
boundary  layers  are  thicker  than  measured. 

Figures  58-61  compare  our  numerical  results  with  measurements  for 
Zwarts’  Mach  boundary  layer.  While  computed  and  measured  shape 
factor  distributions  are  quite  close,  all  other  flow  properties 
show  the  same  general  trend  as  the  other  adverse  pressure  gradient 
cases.  That  is,  the  boundary  layer  tends  to  recover  from  the 
adverse  pressure  gradient  much  more  rapidly  than  measured. 

5.6  FLOW  OVER  TRANSONIC  AIRFOILS 

Continuing  in  our  more  advanced  applications  we  turn  now  to  transonic 
flow  past  two  airfoils,  the  RAE  2822  design  and  the  DSMA  523s. 

For  the  former,  computations  have  been  made  for  five  different  sets 
of  flow  conditions  and  three  different  sets  of  conditions  for  the 
latter.  All  computations  have  been  done  with  our  boundary-layer 
program  EDDYBL  using  measured  pressure  distributions.  To  account 
for  possible  significant  effects  of  streamline  curvature,  the 
Wilcox-Chambers  curvature  term  (Equation  11)  Is  included  for  both 
Surfaces  of  the  airfoil. 

5.6.1  Airfoil  RAE  2822 

For  this  round  of  airfoil  computations,  freestream  Mach  number 
ranges  from  .676  to  *730  while  Reynolds  number  based  on  chord 
length  ranges  from  2.7  to  6.5  million.  All  computations  are 
started  at  the  stagnation  point  and  the  transition  point  adjusted 
to  match  the  measured  location  by  varying  freestream  turbulence 
intensity.  Figures  62—64  compare  computed  and  measured  momentum 
thickness,  shape  factor  and  skin  friction  for  the  upper  surfaces; 


computed  properties  on  the  lower  surface  (dashed  lines)  are  also 
displayed.  As  most  clearly  exhibited  in  the  skin  friction 
distributions,  we  again  find  that  computed  skin  friction  fails 
to  drop  as  low  as  measured  and  the  boundary  layers  all  end  up 
much  farther  from  separation  than  measured.  The  latter  point 
again  suggests  that  the  predicted  boundary  layer  approaches  an 
equilibrium  state  much  more  rapidly  than  measured  as  the  adverse 
gradient  eases. 

5.6.2  Airfoil  DSMA  523s 

Mach  number  and  Reynolds  number  range  from  .6  to  .8  and  from  2  to  4 
million,  respectively,  for  this  round  of  airfoil  computations. 
Figures  65-67  compare  predicted  and  measured  skin  friction  dis¬ 
tributions  for  both  airfoil  surfaces.  On  the  one  hand,  theory 
and  experiment  are  reasonably  close  on  the  upper  surface  for  all 
three  .cases.  On  the  other  hand,  the  measured  boundary  layer  nearly 
separates  on  the  lower  surface  for  each  case  while  the  numerical 
boundary  layer  actually  shows  an  increase  in  skin  friction.  Thus, 
we  again  observe  a  more  rapid  than  measured  return  to  equilibrium 
in  an  adverse  pressure  gradient.  Figures  68-73  display  computed 
and  measured  velocity  profiles.  Except  very  close  to  the  trailing 
edge,  predicted  upper  surface  profiles  differ  from  those  measured 
by  less  than  5% .  Lov/er  surface  profiles  show  larger  differences, 
particularly  in  the  nearly  separated  zone. 

5.7  DIFFUSER  FLOWS 

Thus  far  all  of  our  applications  have  been  to  external  flows  (with 
the  exception  of  the  Samuel-Joubert  case  which  was  treated 
however  as  an  external  flow).  As  part  of  our  overall  objective 
to  cover  as  wide  a  range  of  turbulence  phenomena,  we  now  focus 
on  two  internal  flow  geometries,  viz,  flov;  through,  a  six  degree 
conical  diffuser.  The  first  case  has  low-core  turbulence  while- 
the  second  is  for  high-core  turbulence. 
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5- 7-1  Low-Core  Turbulence 


Figures  74-76  show  computed  and  measured  skin  friction,  velocity 
profiles  and  shear  stress  profiles,  respectively.  As  illustrated 
in  Figure  74,  computed  skin  friction  initially  falls  off  at  about 
the  same  rate  as  measured  but  then,  in  contrast  to  the  measured 
trend,  begins  to  increase  slowly  rather  than  continuing  toward 
separation.  The  velocity  and  shear  stress  profiles  (Figures  75^76) 
show  clearly  that  the  numerical  boundary  layer  ceases  to  grow  as 
rapidly  as  measured  beyond  the  point  where  the  computed  skin  friction 
begins  to  rise.  Again,  numerical  results  suggest  the  numerical 
boundary  layer  heads  toward  an  equilibrium  state  differing  from  the 
measured  state.  In  this  case  it  is  unclear  whether  we  are  ap¬ 
proaching  a  different  equilibrium  state  or  approaching  equilibrium 
more  rapidly  than  measured. 

5.7.2  High-Core  Turbulence 

Computed  results  for  this  case  are  compared  with  corresponding  meas¬ 
urements  in  Figures  77-79-  Although  computed  and  measured  skin 
friction  differences  are  smaller  in  this  case  than  in  the  low-core 
turbulence  case,  computed  skin  friction  variation  again  stands  in 
contrast  to  the  measured  distribution  in  a  similar  manner.  That 
is,  rather  than  decreasing  monotonically,  the  numerical  cf  begins 
to  increase  slowly  as  we  approach  the  outlet.  As  with  the  low- 
core  turbulence  case,  velocity  and  shear  stress  profiles  indicate 
the  numerical  boundary  layer  is  much  thinner  than  measured.  Again 
we  are  either  numerically  approaching  a  different  equilibrium 
state  than  measured  or  approaching  equilibrium  much  more  rapidly. 


5.8  BACKWARD  FACING  STEP 

Our  final  application,  flow  past  a  backward-facing  step,  differs 
from  all  of  our  other  applications  in  a  very  important  way.  Specific 
ly,  this  flow  includes  boundary-layer  separation  while  the  boundary- 
layers  in  all  of  our  other  applications  remain  attached.  Computa¬ 
tionally  there  is  also  a  significant  difference  between  this  case 


and  all  of  our  other  computations,  viz,  we  have  used  ’’surface" 
boundary  conditions  based  on  the  law  of  the  wall  (Equation  15) 
rather  than  integrating  through  the  viscous  sublayer  (Equations 
12-14).  The  use  of  so-called  "wall  functions"  was  necessitated 
by  limited  time  and  funds  for  this  project;  there  is  no  fundamental 
reason,  however,  why  Equations  (12-14.)  can’t  be  used. 

Figures  S 0—83  compare  computed  and  measured  flow  properties  in¬ 
cluding  surface  pressure  distributions,  maximum  shear  stress 
variation,  velocity  profiles  and  shear  stress  profiles*  Predicted 
reattachment  length  is  3*7  step  heights  compared  to  a  measured 
reattachment  length  of  7*0  step  heights.  As  shown  in  Figure  81, 
computed  maximum  shear  stress  is  considerably  higher  than  measured 
through  the  separated  region*  As  the  flow  proceeds  downstream 
of  reattachment  the  numerical  boundary  layer  returns  to  equilibrium 
at  about  the  same  rate  as  measured,  although  substantial  differences 
persist  to  the  final  station  to  which  computation  continues.  In 
this  flow  the  computed  boundary  layer  responds  more  rapidly  to  the 
strong  adverse  pressure  (much  stronger  than  measured)  and  then 
returns  (from  a  quite  different  disturbed  state)  to  equilibrium 
at  about  the  same  rate  as  measured. 


6.  SUMMARY  AMD  CONCLUSION; 


The  model  employed  in  our  computations  predicts  flow  properties  in 
quite  close  agreement  with  experimental  data  for  the  constant-pressure 
boundary-layer,  the  incompressible  mixing  layer  and  for  flows  with 
surface  mass  transfer.  Additionally,  predicted  effects  Mach 
number  and  surface  cooling  on  a  constant-pressure  boundary  layer 
are  close  to  measured  effects. 

For  flows  with  strong  adverse  pressure  gradient,  most  notably  the 
backward- facing  step,  the  model's  predictions  differ  substantially 
from  corresponding  measurements.  As  an  almost  uniform  trend,  the 
model  responds  to  an  adverse  pressure  gradient  In  a  manner  similar 
to  that  which  has  been  experimentally  Observed  Initially,  but, 
upon  removal  of  the  gradient,  returns  to  equilibrium  more  rapidly 
than  measured  (e.g..  Plows  0142,  0143,  8621,  8623). 

Our  success  with  the  flows  with  suction  and  blowing  is  in  part  due 
to  the  careful  research  which  has  gone  into  developing  appropriate 
surface  boundary  conditions  for  such  flows.  This  success  is  perhaps 
an  argument  in  favor  of  Integrating  to  the  surface  (as  opposed  to 
using  wall  functions)  for  this  type  of  flow. 

Our  relative  lack  of  success  in  computing  flows  with  strong  adverse 
pressure  gradient  is  less  easy  to  explain.  Perhaps  we  should  expect 
to  do  poorly  when  the  flow  departs  even  slxghtly  from  equilibrium 


uDon  ooserv: 


;he  gross  discreoaneies  between  computed  and  measured 


Reynolds  stress  development  for  the  homogeneous  cases.  Yet,  this 
would  be  too  easy  an  explanation  as  the  primary  culprit  Is  the  eady- 
viscosity  approximation  in  the  homogeneous  case.  More  plausibly, 
with  a  two-equation  turbulence  model,  we  may  be  attempting  to  des¬ 
cribe  too  much  with  too  little  in  the  turbulent  boundary  layer. 

That  Is*  the  near-wall  portion  of  the  boundary  layer  responds  on 
a  grossly  different  scale  from  that  of  the  defect  layer.  Yet,  we 
attempt,  from  a  single  equation  (the  equation  for  w)  to  deduce  the 
scales  on  which  the  entire  boundary  layer  will  respond  and  change. 


. . 


FIGURES 


All  figures  following  are  in  the  form  submitted  to  the  1981 
Stanford  Olympics.  All  experimental  data  references  can  be 
obtained  from  the  pioceedings  of  the  Conference.  Unless  other 
wise  indicated,  our  computational  results  are  indicated  by 
lines  with  heavy  dots  (  »  •  ) .  Experimental  data  generally 

are  indicated  by  open  symbols. 
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Figure  15.  Homogeneous  shear j  shearing  rate  =  12.9  sec-". 
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Figure  29.  Development  of  an  incompressible  mixing  layer 
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Figure  51.  Mach  2,2  boundary  layer  with  adverse  pressure  gradient. 
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Figure  67.  Airfoil  DSMA  523s;  Mach  .8;  Re  =  3  million 
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APPENDIX 


WALL  FUNCTIONS 

The  purpose  of  this  Appendix  is  to  provide  some  insight  into 
the  use  of  so-called  "wall  functions"  with  advanced  turbulence 
models.  In  its  original  form,  this  Appendix  was  submitted  to 
the  Stanford  Olympics  C  'mittee  as  a  note  and,  except  for  minor 
editorial  changes  needed  -•  or  consistency  with  this  report,  is 
reproduced  in  its  entir-it/. 


1.  Mathematical  Meaning  of  Matching  to  the  Law  of  the  Wall 

Generally  speaking,  in  order  to  solve  the  equations  of  motion  for 
a  viscous  flow  over  a  solid  surface  one  must  specify  boundary  con¬ 
ditions  valid  at  the  surface.  Often,  in  turbulent  flow  computations, 
it  is  convenient  to  avoid  integration  through  the  sublayer.  This 
can  be  done  by  assuming  the  law  of  the  wall  to  be  valid  for  the  flow 
of  interest  so  that  we  write  (for  incompressible  flow) 

1  UTy2 

u2  =  ur(<log  +  B)  (A1) 

where  uT  is  friction  velocity,  <  is  Harman* s  constant,  B~5  for  smooth 
walls  and  v  is  kinematic  viscosity.  The  quantities  u2  and  y2  denote 
tangential  velocity  and  normal  distance  from  the  surface  at  .the  first 
mesh  point  adjacent  to  the  surface. 

The  first  point  I  wishjto  make  is  that  in  a  strict  mathematical 
sense  the  boundary  condition  we  are  actually  using  when  we  invoke 
Equation  (Al)  is: 

i  ury 

u  -*■  uT(ilog  +  B)  as  y  +  0  (A2) 


We  are  in  fact  idealizing  the  flow  as  having  (relative  to  the  over¬ 
all  scale  of  the  boundary  layer)  a  zero  thickness  sublayer. 


The  Origin  of  "Wall  Functions 


When  a  turbulence  model  is  used  which  involves  additional  dif¬ 
ferential  equations  describing  evolution  of  turbulent  field  pro¬ 
perties,  more  boundary  conditions  are  needed.  For  example,  using 
the  Wilcox-Rubesin1  two-equation  model  of  turbulence  we  must  also 
specify  the  appropriate  "boundary"  conditions  for  the  turbulent 
mixing  energy,  e,  and  turbulent  dissipation  rate  per  unit  energy, 
co.  It  is  at  this  point  that  the  concept  of  so-called  "wall  functions 
is  introduced.  These  functions  generally  are  deduced  by  examining 
the  limiting  form  of  the  turbulence-model  equations  as  y+0.  The 
equations  simplify  in  this  limit  primarily  through  dropping  of  the 
convection  terms.  For  example,  in  a  constant  pressure  boundary 
layer  the  Wilcox-Rubesin  model  equations  simplify  in  this  limit 
to  the  following: 


e  du  .  2 

a)  dy  '  ur 


(A3) 


f  <§)  -  e*ue  +  <f  §)  =  0  (Ai 

where  H  is  turbulent  length  scale  defined  by  l  -  e  /w  and  B,  B*, 
y,  d,  o#  are  closure  coefficients  whose  values  are 


B  =  3/20,  B*  =  9/100 

d  =  1/2,  d*  =  1/2 

Y  =  10/9 


(A  6) 


It  is  easy  to  show  that  one  solution  to  Equations  (A3-A5),  which 

we  shall  denote  as  e=e  and  w  =  w  ,  is: 

w  w 


ew  *  ut//gT 
coT.T  =  Mjm  <y 


(AT) 

(A8) 


Qfl 


*- 


where  k  is  the  Karman  constant,  a  fact  which  has  been  arranged  by 
selecting  the  closure  corfficients  to  satisfy  the  condition 
Y  =  $/8*  -  2ok //$*.  Equations  (A7  -  A8)  generally  are  referred 
to  as  wall  functions.  In  computations,  Equations  (A7  -  A8)  are 
used  to  define  e2  and  m2  at  y=y2* 

The  second  point  I  wish  to  make  is  that,  as  with  the  law-of-the- 
wall  velocity  boundary  condition,  the  precise  mathematical  state¬ 
ment  of  the  "wall-function”  boundary  conditions  for  e  and  w  is 


e  *  e 
w  to 


as  y  +  0 


(A9) 


3.  Non-Uniqueness  of  Wall  Functions 

Now,  because  the  equations  for  e  and  w  are  of  second  order.  Equations 
(A7  -  A8)  are  not  the  only  solutions  of  Equations  (A3-A5).  In 
fact,  by  changing  independent  variables  from  y  to  u  one  can  show 
immediately  that  the  e  equation  simplifies  to 


a* 


8*  E2 


1 


(A10) 


2 

where  E  =  e/uT  and  U  *  u/uT-  As  above,  one  solution  has  dE/dU  =  0 
so  that  E  =  1//8*’.  There  are  also  solutions  having  dE/dU  ¥■  0  which 
can  be  obtained  by  multiplying  both  sides  of  Equation  (A10)  by 
dE/dU  and  integrating  twice  to  obtain 


U  =  U 


d  e 


s3  - 


3 

w 


£  +■  A 


(All) 


where  A  is  an  integration  constant  and  EQ,  U0  denote  reference  values 
of  E,  U.  Equation  (All)  is  an  elliptic  integral  whose  properties 
vary  widely  with  the  value  of  the  integration  constant  A.  It  is 
not  my  purpose  here  to  examine  in  detail  the  behavior  of  E  contained 
in  Equation  (All) .  Rather,  I  wish  only  to  emphasise  that  more  than 
one  solution  to  the  model  equations  exists  and  that,  without  care¬ 
ful  analysis,  we  cannot  be  sure  that  all  solutions  necessarily 
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are  consistent  with  the  law  of  the  wall. 

Based  on  this  observation,  the  third  point  I  wish, to  make  is  that 
arbitrarily  deviating  from  Equations  (A7-A9)  by  using  some  value 
other  than  g*  and/or  k  may  introduce  unexpected  surprises,  many 
of  which  may  be  hiding  in  the  integrand  of  Equation  (Ail).  In 
essence,  in  selecting  the  "wall  functions"  defined  in  Equations 
(A7  -  A8)  we  are  (a)  demanding  that  our  boundary  conditions  be 
consistent  with  the  differential  equations  and  (b)  excluding  any 
other  asymptotic  (as  y  0)  behavior  which  might  be  inconsistent 
with  the  law  of  the  wall. 


4.  Effects  of  Pressure  Gradient 

All  of  the  analysis  above  assumes  constant  pressure.  As  will  be 
shown  in  this  section.  Equations  (A7  -  A8)  are  inappropriate 
unless  boundary  conditions  are  applied  much  closer  to  the  surface 
than  is  done  in  common  practice.  To  see  this,  note  that  now 
Equation  (A3)  must  be  replaced  by 


e  du  _ 
w  dy 


U2  +  I|£y 
T  p  dx  ^ 


(A12) 


where  p  is  density  and  dp/dx  is  pressure  gradient.  Equations  (A4) 
and(A5)  remain  as  before.  Letting  y  =  uTy/v  we  can  rewrite 
Equation  (A12)  as 

§§•»  u*  (1  ♦  z&f  y+>  (A13) 

PUT 


Order  of  magnitude  estimates  for  typical  attached  boundary-layer 
flows  in  pressure  gradient  indicate  the  dimensionless  grouping 
multiplying  y+  is  a  small  parameter.  This  suggests  seeking  a 
solution  of  the  form  o 


u. 


e  = 


w  = 


du 

dy 


/gT 

ut 

/Psy 

Uf 


1  +  <j>e^  + 


..) 


(1  + 


...) 


(A14) 


iLi  *  ii.  (1  + 


$u^  +  . . . ) 
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where,  for  simplicity,  we  have  defined  our  small  parameter  <f> 


as  follows: 


_  vdp/dx 


(A15) 


For  the  sake  of  brevity,  I  omit  details  of  the  algebra  and  simply 
state  the  final  solution  up  to  terms  linear  in  <p,  viz. 


36  + 

ei =  fry 

28  + 

W1  =  -31  y 

33  + 
ui  =  -31  y 


(A16) 


Consequently,  for  flows  in  adverse  pressure  gradient  Equation  (A2) 
must  be  replaced  by 


1  ury  qq  *  UTy 

u  -*■  u  (—  log  - — -  +  B  —  - )  as  y  +  0 

T  VK  S  V  31  K  V  .  J 


(A17) 


while,  to  this  order  of  approximation,  the  wall  functions  e,_  and 

W 


uiw  defined  in  Equations  (A7  ^  A8)  must  be  replaced  by 


-•Si 


1  +  36  *  V 
X  31  V 


(A18) 


to  = 

w 


/pxy 


28  .  V 

31  *  v 


CA19) 


To  show  the  importance  of  the  order  <j>  corrections  in  the  wall 
functions.  Figure A1  shows  results  of  a  computation  with  the  Wilcox- 
Rubesin  model  in  which  the  order  $  corrections  were  omitted.  The 
flow  considered  was  the  Bradshaw  "Flow  C"  adverse  pressure  gradient 
boundary  layer?.  In  the  computation  "surface"  boundary  conditions 
were  applied  at  values  of  y+  ranging  between  12  and  20.  The  figure 
compares  computed  turbulent  mixing  energy  at  x  =  7  feet  with 
Equation  (A18).  At  this  position  the  value  of  y+  for  the  mesh  point 


O  COMPUTED  WITH  EQUATIONS 
(A7)  &  (A8) 


EQUATION  (A18) 


nearest  the  surface  is  12.  As  shown,  the  numerical  solution 
rapidly  approaches  the  analytical  solution  given  by  Equation  (Al8) 
and  the  two  solutions  differ  by  less  than  2%  above  y+  =  30. 

In  a  subsequent  solution  in  which  order  <j>  corrections  to  the  law 

of  the  wall  and  the  wall  functions  were  included,  the  numerical 

and  analytical  solutions  are  virtually  identical  up  to  about  y  =70 

2 

beyond  which  point  terms  of  order  <j>  presumably  become  important. 

Skin  friction  for  these  two  numerical  computations  changed  by  less 
than  25?  indicating  the  neglect  of  the  order  <j>  corrections  was  not 
terribly  important  in  this  particular  computation.  However,  please 
note  that  in  applying  the  boundary  conditions  at  y+  =  12  our  error 
in  e  was  only  k%.  Had  we  applied  the  boundary  condition  at  a 
value  of  y  =40  (which  is  typical  for  those  who  use  wall  functions 
in  their  work),  the  error  in  e  increases  to  14?  which  could  very 
easily  result  in  a  skin  friction  error  of  5%  or  more.  Indeed, 

I  did  some  numerical  experimentation  years  ago  and  found  that  cf 
would  change  substantially  with  the  point  of  application  of  the 
boundary  conditions  when  the  order  $  corrections  are  omitted.  The 
fourth  point  I  wish  to  make  is  that,  by  contrast,  solutions  are 
virtually  independent  of  this  point  of  application  when  the  order 
$  corrections  are  Included. 

5.  Other  Effects 

Similar  perturbation  analyses  can  be  used  to  determine  near-wall 

behavior  of  an  advanced  turbulence  model  including  effects  such 

■a 

as  compressibility  and  surface  mass  transfer  .  In  my  work  I 
generally  integrate  through  the  sublayer  so  I  have  found  no  need 
to  derive  any  but  the  leading  order  solutions.  I  have  made  such 
derivations  only  to  investigate  limiting  behavior  of  the  model  in 
order  to  check  for  consistency  with  physical  reality.  Hence,  my 
past  work  offers  no  further  assistance  to  those  wishing  to  use 
wall  functions.  However,  the  procedure  is  no  more  complicated 
than  outlined  above  and  also  in  Reference  3* 
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6 .  Summary 

In  conclusion  I  would  like  to  summarize  the  points  made  above  and 
add  one  further  comment.  The  key  points  I  am  making  are: 

1.  In  matching  to  the  law  of  the  wall  we  are,  in  a 
strict  mathematical  sense,  insisting  upon  specified 
asymptotic  behavior  of  the  velocity/surface-stress 
relationship  in  the  limit  y/6  +  0,  where  6  is  a  length 
characteristic  of  overall  scale  of  the  boundary  layer; 

2.  In  using  "wall  functions"  we  are  likewise  insisting 
upon  specified  asymptotic  behavior  of  turbulence-field 
properties  in  the  limit  y/6  -*■  0; 

3.  Wall  functions  are  not  unique.  For  a  given  turbulence 
model,  there  generally  is  more  than  one  asymptotic 
solution  as  y/6  •+•  0  and  only  one  of  these  solutions 
can  usually  be  said  with  certainty  to  be  consistent 
with  the  law  of  the  wall; 

4.  In  using  wall  functions  for  flows  with  pressure  gradient, 
surface  mass  transfer,  etc.,  solution  accuracy  can  be 
impaired  if  proper  account  is  not  taken  of  these  effects 
upon  the  wall  functions  and  upon  their  point  of 
application. 

If  proper  account  is  taken  of  the  points  above,  it  is  possible 
to  eliminate  at  least  one  key  area  of  uncertainty  in  numerical  work 
and  in  turbulence-model  research  in  which  advanced  turbulence  models 
are  used.  For  example,  as  noted  above,  solutions  become  independent 
of  the  location  of  the  mesh  point  nearest  the  surface  (provided 
y+<60  for  that  point)  when  Equations  (A17-A19)are  used  in  place 
of  Equations (A2)  and(A9)  for  the  Bradshaw  "Flow  C"  case.  Such  an 
end,  I  feel,  more  than  justifies  the  effort  involved  in  performing 
a  straightforward  perturbation  analysis  of  the  asymptotic  behavior 
of  a  given  turbulence  model.  This  behavior  is,  of  course,  unique 
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to  each  model;  Equations  (Al6),for  example,  are  valid  only  for 
the  Wilcox-Rubesin  model.  The  behavior  peculiar  to  any  other 
model  nevertheless  can  be  determined  once  and  for  all  using  the 
same  kind  of  perturbation  analysis.  The  modest  effort  involved 
should  greatly  reduce  uncertainty  about  numerical  precedures 
employed  by  those  who  make  use  of  wall  functions. 
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